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ABSTRACT 


Many variations in design methods for aircraft digital 
flight control have been proposed in the literature. In general, the 
methods fall into two categories; those where the design is done in 
the continuous domain (or s-plane) and those where the design is done 
in the discrete domain (or z-plane). This report evaluates several 
variations of each category and compares them for various flight 
control modes of the Langley TCV Boeing 737 aircraft. Design method 
fidelity is evaluated by examining closed loop root movement end the 
frequency response of the discretely controlled continuous aircraft. 

It was found that all methods provided acceptable performance for 
sample rates greater than 10 cps except the "uncompensated s-plane 
design" method which was acceptable above 20 cps. A design procedure 
based on optimal control methods was proposed that provided the best 
fidelity at very slow sample rates and required no design iterations 
for changing sample rates. 

At very slow sample rates, system roughness or jerkiness 
may be the limiting factor on the design. To better understand this 
phenomenon, an experimental study was defined for the Langley motion 
simulator with candidate designs for evaluation. 
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1.0 Introduction 


Digital Computers for aircraft inner-loop control functions 

1,2 

have been implemented for research purposes and will probably be 
used for operational aircraft within the next few years. This research 
has been carried out in order to facilitate the design and performance 
evaluation of discrete flight control systems for future missions. 

Many variations in design methods exist, but generally fall 
into two broad categories: 1) those where the design is done in the 

continuous domain (or s-plane) and then discretized, and 2) those where 
the design is done in the discrete domain (z- or w-plane). 

Design using the first category is attractive since it utilizes 

the experience gained over many years of continuous autopilot design, 

and the additional discretization step required can be separated from the 

s-plane design process. All the discretization procedures (z-transform, 

3 

bilinear or Tustin's, and z-forms) introduce an error which is small for 
relatively fast sampling (typically > 40 Hz) but grows larger as the 
sample rate is reduced. The effect of this error or approximation in 
the design method is generally checked by a simulation of the system. In 
some cases, a continuous control design exists and the sole requirement 
on the designer is to obtain the software for the digital computer which 
duplicates the continuous control as closely as possible. For these cases 
only the discretization and verification via simulation is required. 

Design using the second category includes the w plane 
4,5 6 

techniques, z-plane Nyquist techniques, and discrete state space 

7 8 9 

techniques. ’ ’ These methods do not introduce an error due to the 
discrete nature of the problem at any sample rate and therefore are 
particularly attractive when attempting to establish the lowest practical 
sample rate. 

There is an interaction between the design method sample rate, 

system sensitivity, roughness of control, prefilter design, and gust 

response that needs to be thoroughly understood in order to synthesize 

9 

digital control systems with confidence. Preliminary results on the 
intci’action between the sample rate and system sensitivity and gust 



response for a high performance aircraft indicate that sample rates 
can be used which are lower than those typically judged necessary in 
previously reported aircraft digital control studies. This 

research primarily concerns design method evaluations as a function of 
sample rate. Preliminary results concerning roughness of control are 
also included. 
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2.0 Design Methods 


2.1 Introduction 

The following sections carry out discrete designs for the NASA TCV 
Boeing 737 autopilot using several different methods. All the resulting 
designs are checked using an exact discrete analysis to determine the 
fidelity of the design method vs. sample rate. 

Although only a few autopilot modes have been studied, the 
dynamics of the aircraft in the remaining modes of the autopilot are not 
sufficiently different to cause any significant change in the results con- 
cerning design methods. However, for other aircraft with significantly 
faster natural modes or greater structural flexibility, the results con- 
tained here do not apply. 

2.2 Digitized S-Plane Design 

Digitized s-plane design refers to the one-for-one replacement 
of continuous transfer function blocks with difference equations that can 
be mechanized on a digital computer. In classical flight control designs, 
these blocks fall into 4 general categories: integrators, washouts, simple 

lag networks, and compound lead or lag shaping networks for frequency 
compensation. The systems being considered here are ones that feed back 
sensor signals (from accelerometers, gyros, altimeters, etc.) which have 
been processed through these networks for stability improvement, or noise 
or bias rejection. 

The digital problem after s-plane synthesis of an appropriate 
control is to obtain the best possible match of a continuous compensator 
with a digital filter. There are many methods of transforming a con- 
tinuous filter into digital form [Refs. 13, 14, 15, 16] among which are 
the direct method, impulse invariant, zero-order hold, bilinear approxi- 
mation (sometimes called Tustin's method), and the matched z-transform. 

The most successful methods exactly map the poles and zeros of the 

st 

continuous filter into the z-plane by the transformatxon z = e for 
all sampling rates. Some methods, notably Tustin's method without pre- 
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warping, cause a shift in the poles and zeros of the filter at low 
sampling rates. Exact mapping of the poles and zeros means that the 
transformation itself plays less of a part in any degradation of 
performance of the digital system. 


The matched z-transform [16] was used in this analysis. For 
the continuous filter transfer function, 

m 


G(s) = 




the digital mechanization is 


G(z) = 


K'(z + l)‘‘.n (2 . e‘“i^)n p cos bjTx + 


rt ^z^ - cos djTz + 


where 

k = r +2S - rn -an, k ^ 0 

and k ' is the normalization constant. k' can be chosen to 
produce the desired gain at either high or low frequency. Thus, to 
choose the gain at low frequency, G(z) is evaluated at z = +1 and 
K picked to yield the desired magnitude. For a high frequency gain, 

G(z) is evaluated at z = -1. Table I shows the transforms and associated 
difference equations for the four forms of compensation mentioned above. 

Even though any root movement due to the transformation is 
eliminated, there is still some degradation of performance due to the 
physical zero-order hold (ZOH). The hold produces a lag proportional to 
the sampling interval, and can be regarded as causing a shift in the z- 
plane zeros of the plant that increases with the sampling interval. This 
zero shift occurs for all plant zeros not at the origin. Thus, any system 
whose digital compensation is determined without taking into account the 
phase shift of the hold mechanization will. have degraded performance for long 
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Matched Z transforms for common filters 



sampling intervals (or low sample rates). 

Formally, the procedure of digitized s-plane design consists 
of transforming an already existing s-plane design by substituting a 
difference equation for each control block. A reasonable measure of the 
success of this technique is to compare the equivalent s-plane root 
locations of the closed loop discrete system at various sampling rates to 
the closed loop roots of the continuous system. This comparison is not 
a function of the merits of the continuous design; the continuous system 
roots are taken to be the desired roots. The consequences of any root 
shifting can be displayed graphically by a transient response simulation. 
This has been done for three systems based on the TCV Boeing 737 : 

1) Pitch SAS, 

2) Pitch attitude control wheel steering, 

3) Longitudinal autoland. 

For each of these systems a block diagram of the model will be presented 
with the difference equations produced by the digitization, as well as 
a locus of root locations. Appropriate transient response plots will 
also be shown. 

The equivalent s-plane roots of the system with a digital 
control law are determined by a method suitable for computer mechanization. 
First, the transition matrix is determined from the plant equations of 
motion and hold characteristics: 

Given 

X = Fx + Gu 

then X , xDx + Fu 
n+1 ^^n n 

^ FT 

where <p = e 

r = \ e^GdT, 

and T = the sampling interval. 

These difference equations are combined with the difference 
equations of the controller (put into appropriate state form), and the 


The models are based on the actual TCV autopilot with nonlinearities 
removed. 
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eigenvalues found. These z-plane eigenvalues are transformed back into 
the s-plane by the inverse relation s = 0n(z)/T. A program to perform 
these operations as well as certain synthesis methods is described in the 
appendix. 

2. 2. 1. Pitch SAS 

The pitch SAS is simple pitch rate feedback to the elevator for 
better damping and higher bandwidth. A short period model of the airplane 
is used here: a block diagram of the continuous system is shown in Figure 1. 

In this system it is assumed that the control law is calculated instanta- 
neously; there is no sampling interval delay at the output of the controller. 
Also note that there is no washout on pitch rate in this model. Although 
a washout is common in most stability augmentation systems, it does not 
affect root locations or their variation since its break frequency is 
quite low. Without the washout, the digital controller has no states and 
no difference equation; therefore, this example is somewhat degenerate 
but will serve as a basis for comparison with other design metnods in 
later sections. The root movement shown in Figure 2 is due entirely to 
the lag induced by the ZOH. This lag causes immediate degradation of the 
system damping. The real root is very slow and is not affected by 
sampling rate. 

2.2,2 Pitch Attitude Control-Wheel Steering 

The aircraft model for the attitude control-wheel steering mode 
is the same as that shown in Figure 1 with the addition of a pitch angle 
state. In this mode pitch rate is fed back through a washout to eliminate 
any steady-state component. The output of the washout and difference 
between the pitch angle and the preselected pitch angle is fed back through 
a one cycle delay to simulate computational delay. 

The existence of this delay depends on the particular mechanization 
and does not always exist. Its presence degrades the stability and per- 
formance of the system more than would be encountered without it. 

A block diagram of this system is shown in Figure 3. 
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The continuous washout for the TCV is: 


d(s) _ 16s 

c(s) ~ 16s + 1 

which yields a discrete transfer function (Table II) of 

d(z) 1 + , z - 1 \ 

c(z) ~ 2 ^ -T/16 ^ 

z - e 

and a difference equation 


d , = C, d„ + C-(c - c ) 
n+1 in 2 n+1 n 


where the values of the coefficients are given in Table II for several 
sample rates. 


elevator short period model 



Short Period model in state space form: 
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Figure 1 Pitch SAS Model 
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Sdesired 


Figure 3 Pitch Attitude Control Wheel Steering Model 


Sampling 

Rate 

Table II Washout Coefficient Values 

-T/16 

C^ = e 

, -T/16 

1+e 

^2 2 

32 cps 

0. 998 

0. 999 

16 

0. 996 

0. 998 

8 

0. 992 

0.996 

4 

0.984 

0.992 


For the purposes of analysis and simulation, the input from the control wheel 
is assumed to be zero. The variation in dominant root locations versus 
sampling rate is shown in Figure 4. The transient response of pitch angle 
and pitch rate is shown in Figure 5 for an initial pitch error of 0. 1 radian. 
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Sampling at 4 cps causes instability and is not shown. 



Figure 4 

Dominant roots of the 
Control-wheel steering mode 

2,2.3 Longitudinal Autoland 

The longitudinal autoland aircraft model is the same as that 
for the control-wheel steering model with the addition of a vertical 
position state given by: 

h = w - 207 . le 

The feedback control (Figure 6) was taken from drawings of the NASA 
Langley's TCW Boeing 737 and linearized. As witb the control-wheel 
steering mode, a one cycle delay is assumed in the controller. The 
digital washout coefficients are the same as in 'liidble II while the 
difference equation for the digital integration is 


^n+1 


= b_ 


^ 2 ^«n>l 


a ) 
n 


- 11 - 


6 (radians) 
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The variation in dominant root locations versus sampling rate is shown 
in Figure 7. Figure 8 shows the error in pitch rate for an initial 
altitude error of 100 feet off the glide slope. Again sampling at 
4 cps results in instability. 

In svimmary, digitized s-plane design is suitable for higher 
sampling rates or for systems with high relative stability. The 
addition of delays due to computation (which are cascaded in several 
multi-computer systems) has an increased effect with slower sampling. 
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2.3 Digitized S-Plane Des ign with Hold Approximations 

A modification to the design method described in section 2.2 is 

to include the effect of the ZOH in the s-plane design. Since this effect 

is sample rate dependent, the s-plane design needs to be refined for any 

selected sample rate. This design technique will be evaluated for the 

Pitch SAS system which was discussed previously in section 2.2.1. In 

this case the s-plane model is represented as shown in Figure 9, where 

G (s) is nominally 1 for the continuous or extremely fast sampling case 
c 

and varied to cause an acceptable design as T increases. 

Assuming that the loss of damping shown to occur in Figure 2 is due 
primarily to the hold characteristics, some compensation such as a 
classical lead: 

H (s) = K (A < B) 

q (s + B(T) ) 

will be required to retain sufficient damping when T increases. The lead 
parameters can be chosen to meet root location or frequency response 
criteria both of which will be described. 


Z G 



Figure 9 Modified s-Plane Design Model 
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2.3.1 Root Location Design 

If the natural frequency and damping of the short period 
roots are preserved, the time response will be essentially unchanged. 
For purposes of constructing s-plane root loci, we approximate the 
transcendental Laplace transform of the zero-order hold by a rational 
expression of polynomials in s. Fade approximates, 

t ™-.m , „vm-l 

a (-sT) + a , (-sT) + — — + a 

m m-1 o 

a (-sT)*^ + a 

n o 

-s T 

of e are formed by matching the largest possible number of terms of 
the respective McLaurin series by suitable choice of the a^'s and b^'s. 
The approximates thus formed are given in Table III (from ref. [17]). 



The success of the approximation depends on selecting a high 
enough order of approximate to match the frequency response of the zero- 
order hold over the bandwidth of interest. Using Table III, we find 
that the exact hold transfer function, 


H (s) 



s 


reduces to: 


and 



Fig. 10 compares the amplitude and phase of all three representations of 
H (s) 
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The ra=n=l approximate is accurate for frequencies less than 1/10 
the sample rate and the m=n=2 approximate is very good for 60 less than 1/3 
the sample rate. For the short period closed loop poles at approximately 
1 cps (Fig. 2), the 2nd order (m,n=2) approximate must be used to evaluate 
sample rates of 5 cps. 

The compensation is designed by arbitrarily locating the lead 
zero at s=-8 and by successive trials determining the gain K and pole 

q 

location -B which restore the short period roots that have been perturbed 
by sampling. This is easily done since the transcendental nature of the 
lag is eliminated by replacing it by a Fade approximate. The results of 
these calculations are contained in Table IV. 

T(sec.) Kg A B 

.05 4.15 8 8.85 

.1 4.1 8 10. 

.2 4.7 8 14.1 

TABLE IV - Compensation Parameters 

Figure 11 shows the equivalent s-plane roots that result after 
converting the s-plane design in Table IV to discrete form using Tustin's 
(or Bilinear) approximation and the matched z-transform. 

2 z*” 1 

The Tustin transform is formed by substituting — — - for s. 

T z+1 

AT/2 - 1 

Kq sj^ Kq (AT/2 + 1) ^ z+ AT/2 + 1 

s+B ^ (BT/2 +1) * BT/2 - 1 

BT/2 + 1 

The matched z-transform was defined in the previous section. As also 
described in the previous section, the equivalent s-plane roots are 
obtained by performing an exact discrete analysis of the system and 
transforming the z-plane closed loop roots back into the s-plane by; 

s = -i ^ z 
T 
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Figure 11 Dominant Root Locations Using Modified 



The root migration in Figure 11 is caused primarily by the 
approximations in the matched z transform and Tustin digitization. 

The figure shows that root migration at the lower sample rates is much 
reduced compared to the identical design (Figure 2) performed without 
the Fade approximate for the hold. 

As a further aid in evaluating these designs, frequency 
respohse data are presented in Figures 12 and 13 for the Tustin digiti- 
zation. The open loop gain = 1 crossover frequency can be seen to be 

lower in the compensated discrete case (GHZ) at to = 5Hz than the 

s 

continuous case (G). Furthermore, the closed loop bandwidth can be 
seen in Figure 13 to decrease as the sample rate decreases. These 
effects occur in spite of the fact that the short period root locations 
remain essentially unchanged; primarily because the additional root 
introduced by the compensation affects the frequency response. 


2.3.2 Frequency Response Design 

An alternate method of choosing the lead compensation is to 
maintain certain characteristics of the open loop frequency response. 
In other words, use the lead to restore the frequency response de- 
gradation due to sampling. 


The closed loop transfer function of plant G with feedback 
compensation H is: 


GZ 

1 + GHZ 

where the zero-order hold (Z) must be included for sampled systems. 


At low frequency the ZOH and lead network have vanishingly 

G 

small effect, and the closed loop transfer function is desired. 

At. high frequency where |ghz| « 1 the transfer function is GZ instead 
of G, and the compensation H cannot counteract the zero order hold. 

But this is well above the bandwidth of the system, the highest frequency 
of interest. The critical location and sharpness of the Band edge, 
(Closed loop magnitude of -3db) is heavily dependent on the frequency 
at which |ghz| = 0 db and the corresponding phase angle. 
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(rad/ sec) 



to (rad/s®c) 

n*. Design Closed Loop Frequency 
Figure 13 S-Plane Design ^ 


s + A 

Thus in the frequency response design, H(s) = K - 

q s + B 

is picked so that the 0 db or cross-over frequency and phase margin 
of G is maintained for GHZ as shown in Figure 14 for 5 cps. 


The lead compensation is designed by arbitrarily locating 
the zero at s = -8sec ^ and choosing K^, and B to meet the phase margin 

and gain criterion of 0 db at = 4. 1 radians/second, the cross-over 
frequency of the continuous plant. The parameters thus obtained are 
contained in Table V. 


Table V 

Lead Compensation by 
Frequency Response Design 


(Sec. ) 

K 

_q 

A 

B 

.01 

4.48 

8 

8.42 

.02 

4.68 

8 

8.88 

.05 

5.41 

8 

10.54 

.1 

7.44 

8 

14.90 

.2 

31.8 

8 

64.53 


In Figure 15 the equivalent s-plane poles of the discrete system 
j using the compensation of Table V are shown. 

The approximate discrete compensation is formed using both 
Tustin (or bilinear) and matched z-transform methods giving rise to two 
different s-plane loci. 

Also shown in the figure are the s-plane root locations before 
digitization . The root movement for the higher order Fade approximates 

is due to the preservation of frequency response characteristics. 

i 

Furthermore, comparison of the various s-plane designs with 
different Fade approximates indicates that the 2nd order Fade yields 
sufficient accuracy for this design although for the higher bandwidth 
roots at T =.2 seconds there is a noticeable difference between the 
2nd and 3rd approximates. 
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Figure 14 Open Loop Frequency Response, 
Frequency Response Design (T = . 2 sec) 
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The discretization changes the predicted s-plane roots 
significantly, especially if the compensation is fast relatiye to the 
sample rate. For our case where A = 8, the Bilinear transform is un- 
acceptably inaccurate at moderate sample rates. The z-transform is 
somewhat more accurate and should be used. The bilinear and matched 
j z-transform and Fade approximate curves coalesce for higher sampling 
and slower compensation. An example of this is Figure 16 which 
resulted from a design using the same frequency response techniques 
but fixing the pole at -9.5 sec ^ and varying the zero for proper compen- 
sation. The compensation is relatively slow, but clearly the curves 
are more consistent in this case. However, if the lead compensation 
pole is too slow, closed loop frequency response will be adversely 
affected. 

; The migration of the short period root that occurs in this 

design (Figure 16) does not adversely affect the closed loop response 
because the resonant peak from the reduction in damping of the short 
j. period is offset by the lag from the real root originating from the 

: compensation pole and zero. Furthermore, the adverse phase lag of 

}1 the zero-order hold is counteracted by the increased bandwidth of the 
! 

I short period poles. 

i The closed loop frequency response characteristics of the 

I designs with lead parameters given in Table V and root locations 

I I given in Figure 15 are shown in Figure 17. Comparing these results 
to the designs based on short period pole fixing ( Figure 13), we 

i 

j find that the magnitude matches better at low frequencies and the 
jj phase is slightly better over all frequencies. 

l;l The designs with slower compensation yielded less discretization 

1 'i 

errors (note correspondence between dotted curves and solid curves in 

i ; 

I Figure 16) because of the slower compensation to be digitized; however, 
the slower compensation yielded a poorer match to the continuous system 
at the lower frequencies than either the fast compensation frequency 
response design on the root location design. 

, j 

iij Either frequency response design criteria or dominant root 

f. location design criteria appear to be feasible methods for adjusting 
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Figure 17 Closed Loop Frequency Response 
For Frequency Domain Design (T = .2 seconds) 
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the s-plane designs for the effect of sampling in the cases studied. 

The frequency response method, however, appears to offer a more reliable 
method in general since movement of all roots is accounted for in the 
design criteria. For flight control modes with higher order dynamics, 
it may be more difficult to pick the dominant poles and adjust the 
compensation without introducing some new dominant pole. 
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2.4 Direct Digital Design 


Direct digital design refers to the design of a controller 
in the z-plane (or the w-plane, with an appropriate transformation ). 

It differs from digitized s-plane design in that the effects of the 
hold circuit are known exactly at all sampling rates. In the z-plane, 
techniques such as root locus and pole placement may be used exactly 
as they are in the s-plane. 

Direct digital design is most widely used in the construction 
of estimators or of controllers with state variable feedback inasmuch 
as optimal control techniques may be applied with little change to the 
discrete system. However, the starting point for most digital systems 
is a continuous mechanization, for which engineers still have the most 
insight as well as an arsenal of design tools. The object, therefore, 
is not to completely redesign the system, but instead, to have changes 
to account for the discrete characteristics and the lag of the hold 
circuit. Toward this aim, two of the systems discussed in the section 
on digitized s-plane design (2.2) were re-analyzed with the addition of 
a lead network at the input to the hold circuit. The lead ratio and 
gain of this network was chosen in order to place the dominant poles 
of the closed-loop system at the location chosen in the continuous 
design. This process can, in general, be done by trial and error; 
but for this effort, a systematic procedure was developed and implemented 
in a computer program (called LEAD, See Appendix A). It calculates 
the lead required from the pole/zero constellation in the z-plane. 

It is important to note that this design was performed in the z-plane, 
which is an exact representation of the digital system. This FORTRAN 
computer program has been expanded so that it contains options to do 
the analysis of the s-plane design method as well as calculating the 
lead network. 

2.4.1 Pitch SAS 

To the SAS system shown in Figure 1 is added a digital lead 
network which performs the equivalent of the continuous lead: 
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using the matched z-transforra method, the digitized version of this 
becomes: 
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The constants for this lead network were determined so as to match the 
equivalent s-plane complex roots of the discrete system with the original 
complex roots at s =-5.15 + 4.5j. Table VI shows the lead network 
characteristics determined and the resulting equivalent s-plane root 
locations. Except for the 4 cps case, the complex roots were matched 
exactly while the position of the additional real root added by the 
network varied. The inability to match precisely the 4 cps complex 

TABLE VI 

Pitch SAS Lead Characteristics 


Sampling 

Rate 



K 


Roots 


32 cps 

8.0 

8.46 

0.917 

-0.77, 

-5.15 + 

4.5, -7.75 

16 

8.0 

9.03 

0.841 

-0.79, 

-5.15 + 

4.5, -7.45 

8 

8.0 

10.7 

0.714 

-0.82, 

-5.15 + 

4.5, -6.75 

4 

8.0 

75. 

0.45 

-0.97 , 

-5.14 + 

3.12,-5.08 


roots was due to the arbitrary design procedure of fixing the lead zero 
at 8 r/sec and varying the lead pole to obtain the desired roots. A 
zero placement closer to the s-plane origin would have allowed matching 
at 4 cps also. 


The additional real root arising from the lead compensation 
dynamics is sufficiently fast to have a relatively small effect on the 
system response. However, for progressively slower sampling the root 
becomes progressively slower and more important in the response. In 
subsequent designs, the detrimental effects of the additional real root 
on the transient response will be evaluated. 
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2.4.2 Pitch Attitude Control-Wheel Steering 


Just as for the SAS described above, an appropriate 
compensation for the control-wheel steering mode is a lead network 
at the output of the digital controller (see Figure 18). LEAD was used 
to match the dominant poles at -1.34 + 2.13J. With decreasing sampling 
rate, however, one of the faster roots moves toward the origin to in- 
fluence the system response. Basically, there is no way to keep all 
the system roots constant with a decreasing sampling rate with only 
one compensation. However, there is a great improvement in performance 
i.e. the response is good down to 4 cps, whereas the system with a 
digitized s-plane design goes unstable. Table Vll presents the para- 
meters of the digital lead network employed, and Figure 19 indicates 
that there is no movement of the root at -1.34 + 2. 13j at all versus 
sampling rate. That figure also shows the position of the next 
faster root at 4 cps versus the lead network zero (sampling faster 
than 4 cps puts this root on the negative real axis and much faster 
than the dominant root). Since these two roots have bandwidths of 
the same order of magnitude, they both must be taken into consider- 
ation in determining performance. Figures 20 and 21 show the pitch 
error and pitch rate transient response for an initial pitch error 
of 0.1 radian with lead compensation (they may be compared directly 
to Figure 5). Although there is some degradation in transient per- 
formance, it is considerably better than that of the system without 
compensation. 



Odesired 


Figure 18 Control Wheel Steering (Discrete Design) 


Table VII 

Control Wheel Steering Lead Characteristics 


imp ling 
rate 

1/T^ 


K 

32 cps 

2.0 

2.23 

0.87 

16 

2.0 

2.49 

0.77 

8 

2.0 

3.31 

0.61 

4 

2.0 

6.03 

0.45 
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Figure 19 Control Wheel Steering Mode Equivalent S-Plane Roots 


Pitch rate (rad/sec) 





2.5 Discrete Design Using Optimal Techniques 

An alternate method to the design techniques discussed so 
far Is to use modern control's "optimal" estimator/controller approach. 

This entails a somewhat arbitrary selection of weighting factors in a 
cost function; however, once these are selected, the discretization 
process is routine. 

An estimator/controller implementation can be exactly the 
same as a classical compensator, the difference being solely in the 
methods used to arrive at the design. In many continuous cases, however, 
the estimator/controller design results in a higher order compensation 
than would result using classical methods; for discrete systems, the 
difference in order is less pronounced. 

Design of compensation using modern methods is made up of 
two parts; the controller and the estimator (sometimes called an observer). 
The controller is determined by minimizing the quadratic performance 
index: 


J 

c 


1 

2 



T T 

X Ax + u Bu 


] 


dt 


and the estimator by minimizing: 


•^o = "i ^ (z - Hx)^ R Nz 



where 

X is the state vector 
X = estimated state vector 
u = control vector 
z = measurement vector 

H = measurement state distribution matrix 
A = state weighting matrix 
B = control weighting matrix 
W = process noise vector 

Q = white process noise power spectral density matrix 
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R = white measurement noise power 

spectral density matrix 

V = z - Hx = measurement error vector 

It can be shown, (ref. [18] [19]) that the solution of the 
optimal control and estimation problem consists of linear state feed- 
back to the control u and linear measurement error feedback to the 
estimator: 

A 

u = Cx 

X = (F + GC)x + K(z - Hx) 

Using eigenvector decomposition methods [19], [20] the 
steady state controller gains, C, and estimator gains, K, can be 
efficiently determined. To accomplish this, the A, B matrices and 
the Q, R values of the noise model are chosen to give a "satisfactory" 
system design. The judgement in determining when the system is satis- 
factory is based on the same factors that one would use for a classically 
designed system. In this particular design exercise, the goal is to 
obtain a modern design which has all the same root locations as the 
previous design methods; thus enabling a direct comparison with the 
other methods. Therefore, A, B, Q, and R were selected to yield 
identical closed loop short period roots while maintaining the additional 
roots which occur with the introduction of the estimator relatively fast. 
The bias estimate state (which performs the same function as the classical 
washout) is an exception to this; its root was maintained at the washout 
break point. 

A continuous controller/estimator for the Pitch SAS with a 
washout was determined for the TCV vehicle. The state space model for 
the aircraft was given in section 2.2.1. Combining this with the 
actuator model and augmenting the state to account for a possible bias 
in the pitch rate measurement (the equivalent of a washout), we obtain 
the system model: 
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where: 


y = Hx + V 
H = [1010]; 

q, = pitch rate measurement bias 
b 

q = pitch rate 

w = z-component of aircraft velocity 

6 = elevator angle in degress 

e 

6 = commanded elevator angle 

e 

c 

y = measurement 

The values of A and B selected for use in the cost function 


were: 


A = 


0 

.33 

0 

0 


0 

0 

7 . 3sX^ 
0 


0 
0 
0 

-.86 


B = 1 


which yielded: 


C = [0, -1.3^‘®, 247.23, 9.(N"® ] 
the values of Q and R selected are: 

r 

Q = 


which yielded: 


22.877 

0. 

0. 

0. 

0. 

5.0^® 

0. 

0. 

0. 

0. 

33.92 

0. 

0. 

0. . 

0." 

\7 

1.106^ 

V. 

r 

1 



R = 1 


K = 


4.783 

-i.sos'^ 

8.816 

-314.5 


r 

the notation 3.1> 


means 3.1 x 10 


-3 
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For comparison purposes, the poles and zeros of the resulting 
design’s open loop transfer function were computed . Figure ( 21 ) shows 
the selected closed loop root locations and a root locus vs the open 
loop gain. This design yields three roots more than a classical design 
with a washout. If compensating the system using classical methods to 
overcome effects of sampling, one or two extra roots may be added, a 
step that is not required when discretizing a controller/estimator design 

By formulating the system equations in a discrete fashion but 
using an identical cost function and weighting values to the continuous 
design, an optimal discrete controller/estimator results which closely 
matches the performance of the continuous design. 


In this case, the discrete system representation is as stated 
in section 2,2 and results in a discrete controller: 


+ K(Z. -Ifa. 1 
1+1 1+1 


or 


u . 
1 

= Cx. 
1 

A 

^i+1 

= X. , 
1+1 

^i+1 

■i- 

II 

z . 
1 

= Hx.+v 
1 

X . , 

1+1 

II 


A , 

- X. , 
1 1 


i+1 


The actual cost function used in the minimization was discrete, 


1 . e. 


00 ^ _ 
1 \ P* T X 

Controller: min J = — S lx. Ax. + u. B u. I 

u. 2.„Lx Di 1 DiJ 


i=0 


Estimator: 



A T 

Hx . ) R (z . 
1 D 1 



however, the values of A and B are computed based on A and B 

D D 

from the continuous design. The values of and are the same as 

the Q and R used in the continuous system optimization. Eigenvector 

[19,201 

decomposition methods yield a very efficient tool in solving 

for the optimal gains for either the continuous or discrete problem. 
Figure (22) shows the equivalent s-plane root locations for the discrete 
system designed as a continuous system in Figure (21), The roots show 
very little migration due to sampling. Figure (23) shows the closed loop 
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Figure 22 Root Locus of Controller/Estimator Design vs. Sample Period 



to (rad/sec) 



Figure 23 r Closed loop Frequency Response versus 
Sample Period for Design Using Quadratic Synthesis 




frequency response characteristics for the various sample rates. To 
demonstrate the modification to the compensator that is being brought 
about by this digitization procedure, we have included Figure (24). 
Note the very large changes that occur in some gains vs. sample rate. 

The digitization procedure is essentially routine and no 
design iteration is required for any of the sample rates evaluated. 

The design iteration is carried out on the continuous system with 
the end product being a set of controller weighting matrices, A & 

B, and estimator noise matrices Q & R. ' These matrices are easily 
turned into discrete control laws (or compensators) for any reasonable 
sample period in a non-iterative manner. 

The most difficult part of the design method is obtaining 
the continuous optimal design. For systems being designed for the 
first time, the iterations required are on the same order as those 
required for a classical design. For systems whose continuous design 
exists and was obtained using classical means, the iterations required 
to obtain a matching optimal design is somewhat cumbersome, but the 
requirement to "match" only exists in a study such as this where it 
is important to make consistent comparisons. 


actually can’t be matched exactly because the order of the continuous 
compensation is higher for the optimal method. 
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Controller Gains 




Figure 24 Controller and Estimator Gains Versus 
Sample Period 



2 . 6 Comparisons of Design Methods 


Five different design methods have been evaluated for various 
modes of the TCV flight control system. For a common comparison of 
the methods, the fastest loop (inner pitch stabilization) has been 
evaluated by all of the methods. 

Figure 25 shows the equivalent s-plane short period root 
locations of the discrete systems for four of the design methods 
and Figures 26 and 27 show the closed loop frequency response character 
istics for 5 of the methods. The methods are: 

1) Uncompensated s-plane design (section 2,2) 

2) Modified s-plane design based on pole fixing 
with the matched z-transform digitization 
(section 2.3.1) 

3) Modified s-plane design based on frequency response 
with the matched z-transform digitization(section 2.3,2) 

4) Discrete design based on pole fixing (section 2.4) 

5) Optimal discrete design (section 2.5) 

Figure 25 shows the most severe degradation in stability for the 
"uncompensated" design method (as expected), a conclusion that is 
amplified by Figures 26 and 27. For this flight control design 
problem, the uncompensated design method should not be used with sample 
rates much slower than 20 cps. Since the system dominant poles are 
at approximately 1 cps, this represents sampling at 20 times the 
dominant closed loop roots. 

Further comparison of the designs shows a large root movement 
of the frequency response design but no corresponding change in the 
frequency response. This arises because the short period root movement 
is being counteracted by shifts in the other real roots which are not 
shown and demonstrates the usefulness of the frequency response 
information in performing a total design evaluation. 

With the exception of the uncompensated design, all other 
design methods would probably be considered acceptable; however, most 
show some loss of bandwidth at the slow sampling rate (5 cps). 
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Figure 25 Equivalent S-Plane Short Period Root Locations 
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Figure 26 Comparison of Closed Loop Magnitude Frequency 
Response by Design for T = .2 Sec. 








Since the frequency response design was based on retaining 
the same gain and phase at crossover ( 4.1 r/sec) it’s interesting to 
note that the closed loop frequency response information does not 
match the continuous system at 4.1 r/sec. This is due to the error- 
brought about by the digitization procedure (the matched z-transform) , 
an error that is present for any design method that is based on 
s-plane results, regardless of whether the hold characteristics are 
accounted for or not. 

Note also however, that accounting correctly for the 
digitization process is not nearly sufficient in itself to obtain a 
good design. This is demonstrated by the discrete pole fixing design 
in Figures 26 and 27, which, except for the uncompensated design, is 
the worst match to the continuous of all the designs. This design was 
not shown on the root locus in Figure 25 because it was based on 
maintaining the z-plane roots at the location which produced no move- 
ment of the equivalent s-plane roots. The shift in the frequency 
response characteristics is due to shifting real root locations which 
were not taken into account when compensating the system. 

Factors other than the accuracy of the method must also be 
taken into account when evaluating design procedures. Table VIII 
is included to illustrate these differences. The first column refers 
to the quality of match to the continuous closed loop bandwidth as 
evidenced by Figure 26. Intuitive design refers primarily to whether 
the design iterations are performed in the s-plane or not, which most 
designers prefer; however, the modern discrete is given a qualified 
yes because many designers do not consider the weighting matrix 
selection an intuitive process. The third column refers to whether 
the digitization effects are evaluated in the z-plane or not. For the 
s-plane methods, the magnitude of the error is dependent on the 
frequency placement of the compensation — the higher the frequencies, 
the worse the error. This error can always be computed by an 
analysis such as that performed by this study, but this is not con- 
sidered to be part of the design method. The fourth column refers to 
whether any design iterations are required when changing sample rate. 
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a feature that is attractive if part of the design exercise is to determine 
the minimum practical sample rate; hence designs at many sample rates being 
required. The order of the compensation refers to the complexity of the 
control mechanization. For all but the modern design, considerable 
freedom exists in the order of the compensation, although the order 
tends to increase with slower sampling. For the modern design the 
compensation order is typically equal to the plant order but can be made 
equal to the plant minus the number of measurements when using a reduced 
order estimator. The last column is a somewhat arbitrary attempt to 
indicate at what sample rates the method can be used. For example the 
uncompensated design method should not be used at sample rates below 
20 times the dominant closed loop poles. These estimates are based 
primarily on the quality of the closed loop bandwidth match as indicated 
in Figure 26 and in the more detailed design data in the individual 
sections. They are valid for the example used in the designs, the 
737 TCV pitch stabilization mode, and should be used with some caution 
for other design problems. 
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3 . 0 Roughness 


Since the output of a ZOH is constant for an entire sampling 
period, the control is a series of steps which may be perceived by a 
pilot through his motion or visual cues as jerkiness or roughness. 

These effects will be most noticeable in the high frequency modes of 
the aircraft and become increasingly noticeable as the sample rate 
decreases. 

To evaluate and quantify this effect of digital flight 
control for varying sample rates, the following designs have been 
selected for implementation on the Langley motion simulator. The 
autoland and control wheel steering modes have been selected for this 
investigation. The autoland mode essentially places the pilot in a 
monitoring type of function and may produce different reactions compared 
to the control wheel steering where he will be required to execute 
various tasks. 

The holds to be considered are (Fig. 28) the standard zero 
order hold (ZOH) for a baseline and the integral hold. 



Figure 28 Hold Definition 
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The integral hold (IGH) is motivated by the desire to 
eliminate the step changes in u characteristic of the ZOH but to 
retain its hardware simplicity. Its mechanization simply requires 
the addition of an anolog integrator between the standard ZOH and 
the actuators in the aircraft. 

Figure 29 combines the different elements of the evaluation 
and indicates where the different portions of the computations are to 
be carried out . 

The algorithms to be implemented for the hold are: 


ZOH: 

u(t) 

= u. ; t . < 


1 1 




IGH: 

u(t) 



The aircraft model and autopilot for the control wheel 
steering mode for use with the ZOH is given in Figure 18 with the 
washout coefficients given in Table II and lead coefficients in 
Table VII for the sample rates to be used on the CDC 6600, 


The aircraft model for the autoland mode is given by 
Figures 1 and 3 and the additional position state relationship: 
h = w - 207 .1 0 


The autoland autopilot for use with a ZOH is given in Figure 6. The 
digital washout coefficients are again given in Table II and the 


integration by: 


+ 2 


For adequate stability with decreasing sample rates, the overall loop 
gain shown to be 0.85 in Figure 6 was adjusted according to Table IX. 


TABLE IX 
Autoland Gains 


Sampling Rate Gain 

32 cps 0.60 

16 0.51 

8 0.34 

4 0.21 
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Analog Computers (EAI 231R*s) 



Digital Computer (CDC 6600) 


Figure 29 Roughness Simulation Schematic 
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For the integral hold, designs are preliminary and need 
further refinement for good response characteristics. Because of 
the extra 90° phase lag of the Integrator and the additional feed- 
back from the integrator, the design is more difficult than that 
for the ZOH. One preliminary design was performed using a modified 
s-plane approach where an approximation to the delay is contained 
in the system model. The design was obtained for 6 feedback only 
so did not represent the full control wheel steering mode. Figure 30 
shows the block diagram and defines portions of the compensation. 


I 



digital ^ I ^ analog 


Figure 30 Integral Hold Design 

The transient response of this controller obtained from a 
hybrid simulation with a 10 cps sample rate in the digital portion 
is shown in Figure 31. It shows that the performance can be made 
to match that of the system without the extra integrator even though 
the controller input to the actuator is ramped instead of stepped. 
Furthermore, the magnitude closed loop frequency response is shown 
in Figure 32 to match very closely as well. Further design efforts are 
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Figure 31 Transient Response of Pitch Stabilization Mode 

from Hybrid Simulation 

underway to determine compensation for the control wheel steering and 
autoland modes. 
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Figure 32 Closed loop Frequency Response Comparison 
of the Integral and Conventional Holds 

Analytical evaluations of roughness have been proposed [20] 
and can be used for a quantitative comparison with pilot evaluations; 
they are: 


where the terms are defined in Figure 33. These functions are to be 
evaluated in the digital computer during the simulations. 

The exact nature of the evaluations with regards to initial 
conditions for the autoland, test maneuvers for the control wheel 
steering, and form of pilot evaluations, have not been determined. 


ZERO ORDER HOLD 


Yi -^2 

T 2T ^ t 




Figure 33 Roughness Function Parameter Definitions 
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CONCLUSIONS 


Discrete Design methods have been compared at various sample 
rates for several different flight control modes for the Langley TCV 
Boeing 737. It has been shown that digitization of continuous auto- 
pilot design leads to a substantial degradation in performance for 
sample rates slower than 20 Hz (or 20 times the closed loop bandwidth). 
Modifications to the s-plane design to account for the lag due to 
sampling improves performance considerably and can be used with 
confidence at sample rates of 10 Hz (or 10 times the closed loop 
bandwidth). For sample rates much lower than 10 Hz, a discrete design 
procedure is advised to eliminate the digitization error present in 
any s-plane design procedure. A new discrete design procedure has 
been evaluated which is based on optimal control and estimation 
techniques. It entails performing the design iterations on the 
continuous system, then using the resulting weighting matrices in a 
discrete optimal control procedure. The resulting discrete designs 
match the continuously designed system better than any other method 
evaluated and requires no iteration during the discretization 
procedure, a feature that no other acceptably accurate method possesses. 

A test of roughness has been defined and control algorithms 
determined for the Langley simulator facility with sample rates of 
32, 16, 8 & 4 Hz, Preliminary designs are Included for the integral 
hold, a hold which requires the addition of an integrator between 
the ZOH and the actuator. 
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Appendix A: Program LEAD 


Program LEAD calculates the eigenvalues in the z-plane of a 
continuous system with a ZOH, actuator dynamics, and a digital controller 
mechanized in difference equations. An option of the program allows the 
user to calculate the lead network required at the input of the ZOH to 
produce given s-plane roots. 


Derivation; 

A system may be written 
X = A X + B u 


A1 


where 



a state vector partitioned into the original 

states X , from x = Ax + Bu, and x , the 
s a 

actuator states 




A2 


Thus u is the control signal to the actuators whose differential 
equations are of the form 

X = Ax + Bx A3 

s a 

If the transition matrix is calculated for equation A1 


where 


X = PX + Qu 
n+1 n n 


P = e 


a't 


A4 


Q = 


A V / 
e B dT 


"^o 

To equation A4 must be added the difference equations of the digital 
controller. The control is of the form 


u 

n 


= C,x 
1 n 


C„x 

2 cn 


+ C„x 
3 sn 


A5 
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where 


X is the controller state variables, 
c 


then 


X is the derivative of original state variables such as 
s 

acceleration, (and x = A x) 

s 


X - = [P + Q(C + C A )] X + QC X 
n+1 1 3 n 2cn 


A6 


The difference equations can be written 

• • ' 

X , = M X + M„x , + M„x + M,x , + M^x^ + M„x , A7 
cn+1 1 cn 2 cn+1 3 n 4 n+1 5 5n 6 sn+1 

The presence of x , , x , , x , is due to the existence of these 
cn+1’ n+l’ sn+1 

types of terms in most matched 3-transforms. 


Now 


X , = A X 
sn+1 n+1 


which can be rewritten using A6 
and finally 


cn 


X , = (I-M„)"^ j M„ + (M^ + M a')[P + Q(C, + C-A')] + M. A ' > x 
cn+1 2 ^3 4 6 13 5\n 


+ (I-M^) 


■‘[“i 


+ (M, + M^A 
4 6 




cn 


A8 


A9 


AlO 


Thus, given the state definition A2, A3, the sampling interval T, the 
control law A5 , and the difference equations A7 , LEAD solves for the 
eigenvalues of the matrix 


p + q(Cj^ + CjA ) 


<JC„ 


♦ (M^ + Mg A')tP + Q(C^ + QCj) + + (M^ + M^a' )QC j] 


To calculate the lead network required to fix a specific root 
in the s-plane, the program first computes the poles and zeroes of the 
system in the z-plane. It then transfonns the desired root into the 
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z-plane, calculate the lead angle deficiency, and places the pole of 
the lead network from this. The gain is computed by standard root- 
locus methods. 

The user must program 3 FORTRAN subroutines: size, diffeq, 
and plant. They define respectively the size of various arrays, the 
difference equation matrices, and the state equations of the plant 
and actuators. A listing follows of the main program and relevant 
subroutines. 
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IMPLICIT REAL (K,L,MJ 
REAL^-3 F, 1CNT,INT,GN,C9,ACL 
REAL*8 SRE,SIM,ChI,CWR 
RfAL*6 CLCGfOSORTfCATAN 

01 MEN SICA ICmT (20) , 1NT(20 ) , GN ( lO ,1 j ) , 09 1 ^O ) » Cv. M lu) , C WI ( 10) 

5 fTIUVUO.lO) ,ACL(1C,1'0) 

01 MS A SICA RVtC (10,10) , \/ FCRA ( 10 , i 0 ) , VECiNdO, lu ) , VECR UO , 10) 

1 , VECI ( 10,10) ,VRV( 10) ,VI V (10) , VHRV ( lO) , V/i 1 V 1 1 0) 

OIMEhSICN APRIME(10,10I,OPPIMEUu,iu) ,P( iO,iU) ,g(iu,10) 

1 , Cl (10,10) ,C2< 10,10 ) ,C3 (10,10) , Mi (lO, 10) ,M2( 10,10) 

2 , M3 ( 10,10) ,M9( 10,10) ,M5 ( 10,10) ,M£> ( 10,10) 

3 ,A( 1G,10) ,B(10, 10 ) ,C(10, 10),0(10,lO) ,0.1 uMU 0, lu) , TcM.3BA( 10,10) 

4 , TEMAA(10,10) , TEM2AM10 ,1C) ,Tfct-lBA ( lO , 10 ) , T uM Ad ( 10 , 10 ) 

2 , VI (20 ) ,V2(2C) ,TEM3EA( 1C ,10) ,TtMJa( 10,10), Ti:,-I^;BB( 10, 10) 

6 , TEM2AE( 10,10 ) ,F( 10,10 ) ,TEH2BA( 10,10) , TcM4BA ( lO, 10) 

OAT A APRIME/U)0«Q. / ,BPRIM£ / 100*0./, P/1 J0*0./ ,w/iOO*J./ 

3 , Cl /10C*0./, C2/10C»C. /, C3 /100*0. / ,Mi/i00*0. / ,.M2/i0U*0./ 

2 ,M3/100»0./,M4/100»Q./,M5/100*0./,M6/iOO*0./ 

CAT A A/1CC*0./,B/1CO*0./ ,C/ 100*0./ 

CALL SIZE(APLAKT,A0IG,AC,AACT) 

NS=APLANT4-NACT 
NSTCT=NS+NDIG 
NST0T2 = 2*.NST0T 

CALL CniT(APLANT,N0IG,NC,NACT,Nd,ilSTQT,UST0T2, 

1 ICNT,IHT,GN,09,CV,R,CU,T lAV, ACL,RVcC,VEv.rti<, VtClN, VECR, 

2 VECI ,VRV, VIV, VRP.V, VR IV.AFR IMf ,6PtU Mc,P,Q,Ci,C2,C3,Ml ,M2,M3,M4, 

3 M5 ,M6, A,e,C,D,TF.M58A, TEMAA,T£M2AA, Ti.3iiA , TtHAb ,Vi,V2, 

4 TEM3eA,TEMaB,T£M23B,16M2Ae ,F, TEMZBA, TEM4BA, D.IUM) 

RETURN 
END 

i SLESGUTINE 001 T ( NPLANT , NO 1 G ,NC ,I1ACT ,dS , NSTOT ST0T2 , 

1 ICNT,INT,GN,09,CV»R,CV<I,TINV,ACL,RVtC,V£CR.<,Vi:CIiJ,VtCR, 

2 VECI ,VF V,VIV,VRRV,VRIV,APRIME,bPKiM2,P,W,Cl,L2,C3,Ml,M2,M3,M4, 

i. 3 M5,ME,A,e,C,0,TEM5CA,TEMAA»T£M2AA, iEMBA,r£MAB,Vl,V2, 

4 TrM3EA,TrMBB,TEM298,TtM2Ae,F,TEM2BA,TcM43A,0liUM) 

IMPLICIT PEAL (K,L,M) 

•REAL* 8 F, ICNT,IMT,GN,CS,ACL 
P.r,AL*8 SRc,SIM ,CkI ,CWR 
REAL*8 DLCG,OSOkT,OATAN 
CCMM.cn /T /T SAMP, CPGLES ,C/E PCS 
CIMENSICN ZtRa(2,8) ,PCLE( 2,8) 

OIMENSICN ICNT(NSTGT2) ,I.VT{MSTCT2),Gta\jT0r,Hi UT) ,09(NST0T2), 

1 CViRiNSTGT) ,CWI (M5TGT ) , TINV (MC IG, NJ 1 0 ) , ACl lU SI OT , NS T J T ) , 

2 RVEC (NSTCT,f)STG T ) , V t CT N ( ^ S TOT ,,)S7 j T ) , VcC In ( hS TuT , NSTCT ) , 

3 VECR (NSTOT.NSTCr ) , VEC I ( NS TCT ,.NST jT ) , vR V( i,o I CD , VI V ( ,'( STOT I , 

4 VRRV (.\STCn , VRI VINSTCT ) , AFSl Wc (iiS , ilS ) , dPRlMEI nS , NC ) , 

5 P(r.S,MS),0(.'Ji,NC),Cl(iJC,NS),C2(UC,i,j;G),C3(iiC,:ij) 
oiNcr.'srcN Mi{r:r!!r,,'!nrr ) ,Wi(NoiG,r; ji Ot-idCuoiG.NS) ,m4{noig,ns? , 

1 H5(NC:G,nS) ,M6(,‘)CIG,NS ), A ( NPLA;iT , GPlAnI ) , B( riPLANT ,NC ) , 

2 CLMACT ,NAC r ) ,D(NACT ,NC ) , TEM5P,A(NDiG,.lS) . l£:rtAA (.'JS,;(S) , 

3 TcM2AA(NS,MSl,TEMRA(NCIG,NS),TEHAti(i',S, flD IG) (ViInSTUTZI , 

4 V2 l.•>.STl:T2) ,TEM3PA(;iC IG,.N'S ) ,TSM3B( .,0 1 G, ,-,0i G ) ,T c.-M2dB ( .■^ U I C, NO IG), 
3 rL?''2AB (.NS .NDIG) , F { r;S TC T , N £ TU T ) ,T;:M 23A( ilul mS ) , Tc ( MOIG ,NS ) 

6 ,ONC M( i.ACT ,NC ) 

CALL PLANT(NPLANT,NOIG,NC ,NACT,fiS,i’4Sror,i,S rOTc, 

:i !C4T , IGT ,CV,:;,Ck I ,T INV, AGL,,\V tc, VECrSN, VcClN, vECi^, 

2 V'J.: I , V;. V , Vi '/. Vr-,RV , ./•■■. IV ,APr IME,ul'.Ai M;; ,P ,(j,Cl ,C2 ,Cj, .U ,M2,.*.<3,M4, 

3 '• V : M t , A , fi , C 1 i3 , 1 1 : M 3 t: A » f C >? r , T t . I 3 NBA , T , V i , V ,1 , 

4 1 i.. M .5 L -A , 7 t M [j ( j X T 1; M 2 B i; , f , F , T ^ - .2 0 t u .M4 A f , ..4^1 ) 

IGU'? Ci\LI. C'PCMFT ( ' ZcRG? ’ , £ , ).C ) 
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99 




190 J 
1909 


REAC CCNFENSATION ZERC (PCSITIVt SlGNl 

REACdO,*) CZEROS 

CCNTINUE 

CALL FRC^'PT( • INS? 

READ INSIPUCTIONS 

INS = l,2f3t OF. 9 (INTEGER) 

TSAMP = SAMPLING INTERVAL (SECONDS) 

TEE THIRD AND FOURTH VARIABLES READ HAVE rtcANlNGS DEPENDENT 
GN THE VALUE OF INS 

INS=1; RCCTS IN THE Z PLANE AND S PLANc WILL SE FOUND FOR 
TFE SYSTEM MTh AN INPUnED GAIN A.,0 COMPENSATION 
PCLE AMD ZERO. THE THIRD VAKImLUE IS Trie GAIN. IT 
.MULTIPLIES ALL THE ELEMENTS uF THE UNuM HAIR IX. THE 
FCURTH VARIABLE IS THE CuMPEMSAT I uN PULc (PDSiTIVE 
SIGN). IF NC CCMPENSATION IS DESiRED, SET THE POLE 
LOCATION AT THE ZEFC LCCATIOli. Trie RODTi WILL THEN 
INCLUDE AN UNCCUPLFC ONE CDRReSPOUUlNG TD THE 
CCMPENSATICN. 

INS=2; THE GAIN AND PCLE LCCATIQN WILL Be CAlCULAFcD THAT 

WILL GUARANTEE PARTICULAR ROUTS IN THL S PLANE. THE 
DESIRED S FLAME RGCT 15 luPUTTEO AS THc THiRO AMD FOURTH 
VARIABLES (REAL AND I.MAGINARY). AS A CHcCK, THE 
fiOCTS APE CALCLLATED WITH THIS COMPUTED CDmPE.NSATIOM AND 
GAIN. THIS CAN BE LSED TO CHECK. UN Trie LOCATIDN OF 
CTHER roots WhCSE LCCATIuNS WERE nuT SPECIFiED. 

THIS OPTION WILL WCRK CURRENTLY ONLY FDR 1 ACTUATOR 
THAT INFLUENCES THE STATES THROUGH Ohl nJW-ZERU ELFMEMT 
CF DNQ.M. ANY CTfER FORMULATION WILl Gi Ve UNPREDICTABLE 
RESULTS IN THIS VERSION. 

INS=3; this performs the same CuMPUTATION mS FuR InS=2. 

THIS CAN BE USED ORLV AFTE.R INi = 2 HAS BEEN SPECIFIED 
FCF A PARTICULAR ZERO LUCATIOn AND iAMPLiriG INTERVAL. 
THIS OPTICM EXISTS TC ELIMINATE SO.NE CCMPUTATIOM THAT 
WAS PERFORMED FREVICUSLY. 

IMS = 9; RETURNS TO A .‘’PCMPT FUR A NLK ZcRU. The SECOND, THI.RD 
AMD FOURTH VARIABLES MUjT BE INPUTTcD, BuT ARE IGNORED. 


READdO,*) INS,TSAMF,SfiPN7,SIPNT 

GO TO (99,99,166,1009 ), INS 

CGA IN=SRFRT 

CPOLES=SIFNT 

GO TO 199 

GAIN=1. 


CPGLES=CZEROS 

CALL CMAT(D,0N0M.,«ACT,NC,GAIN) 

CALL Cr;TP.L( NPLANT,r.!DIG,MC ,NACT,NS,nSTJT,.NSTDT2 , 

ICNT , INT ,GM,D9,CWPv ,CW I ,T lAV, ACL ,RVEC , VECRn, Vr.CIN, VECR, 

V I; C I , V .-. f V I V , V i '. f. V I V . I V , A r R 1 1-1 E , i> P .\ I F. e , P , 0 f V X , e , , ri jt , A , M 3 

M3,M6 ,A ,B,C,D, FEV5BA ,TEMAA ,TEM2AA,TcMBA,TLMAt3 ,Vi,V2, 

TE M3BA, TEMBB, TE»'2 3B, TFM2 AB ,F ,TtMEBA, leMAoA, DMuM) 



CO 1903 J=1,NST0T 
P = Cwj-. (J)^»2*^CWI( J)’'*2 


IF (R .GT. 1. ) GC TG 1903 
XR=CWR( J ) 

YI=CWI(J> 

GO TO 1909 


CCMTIM.E 


GAIM=C. 

CALL C‘-'AT(D,DNCM , MAC r,isC, GAIN) 

CALL CNif'.L (IMPLANT ,ND I G , NC , f- AC T , NS , . JSTo f , <«S TUT2 , 
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1 ICNT,INT,GN,n9,CWP ,CM,TINV,ACL,RVLC,VcCkNi VcCINf VcCRt 

2 VECI ,VRV,VlV,VFRV,VRIV,APRrMF,BP^lMc,P»U,i.i,L2,C3f Mt ,N2,M3,M^, 

3 P5 ,P.6 , A. e tC,0, T£K5BA ,TEHAA,TQ.-iiAAt TiiMBA , Tti'lAB ,Vi, V2» 

. A TfcM3BA,TEMBB,TRM2BB,7EP2AE,F,TEM2bA,TtM43A,CriuM) 

SF=1. 

CO AOA I=1.NSTQT 

SF=5F/0SCPT((XP-CWR( I > )**2*(YI-CW1 ( in*»2) 

PCLFd, I )=CWR( I) 

AOA P 0 LE( 2 , r )=cwi { n 

GAIN=1.R10 

CALL CMAT(0,DhOP,NACT,NC,GAIN) 

CALL CNTF,L(APLANT,NCIG,NCtNACT,Ni,Nl>TL)T,NSTUT2, 

1 ICNT, IM,GN,09, CRR , CU I , T IN V, ACl ,K VtC t ViCRu, Vt CiN, \/i;CR. 

2 V5CI ,VRV,VIV,VRP.V,VR I V , A FR I HE , BPAI Me ,P , C , Ci , C 2» C3 , Ml , M2, M3 ,MA, 

3 M5,H6, A,a,C,D,Tew5BA,T£MAA,TFM2AA, TeMoA, TcMAd ,Yi,V2, 

4 TEM3EA,T£MBB,TEH2eB,TEP2AE ,F,TEM23A,TEM4dA,D.JLiM> 

NZ=0 

CO 505 I=1,NSTCT 

IF (Ck^Rd »>^‘*=2*CWI( n*^2 .GT. 2J GO TJ 505 
NZ=NZ+1 

ZERC( l.MZ )=CWRd I 
Z£R0(2,NZ)=CWId ) 

505 CCKTINLE 

CG 609 1=1, NZ / 

609 SF= <F*SCRT( {XR-ZefiOd , I » YI-Z£K0(2, I n*»2 » 

166 C ZeROZ=eXP {-rSAMP>^C2EP0S » 

ZPPNT = EXP( SRPNT’i'TSAMP )-»CCS(SIPnT*TSAMP) 

ZlPNT = 6XP(SRPNT*TSAMPI*SIMSIPNT»foAMP) 

CALL LCACIZRPNT, Z1PNT,NSTCT, POLE, NZ, ZERO, CbAlM,CPQLtZ, PHI 
1 ,CZEPJZ,SF) 

CPOLES=-ALOG(CPOLEZ) /TSAHP 
WRITE (6,£440> CPQLHS.CGAIN 
5440 FOPMATd POLE •,G15.6,* GAIN *,G15.6) 

199 CCMTINUe 

CALL CMAT(0,ONGM,NACT,NC,CGAIN) 

CALL CNTRU NPLAriT,N'OIG,NC ,N ACT , Nj , NoTuT .MoTOfZ , 

1 ICNT ,INT,GN,D9,CNR,CW1 ,TINV,ACL,RVEC,VECRN, VdCIN, VCCR, 

2 VEC! ,Vf.V,VIV,VHRV,VRIV,APPIME,0Pr<IHE,P,O,Ll,Cz,C.3,Mi,M2,M3,M4, 

3 V.5,M6,A,e,C,0,TEM5B<V,TFMAA,TCM2AA,TcMaA,TcMAd,Vi,V2, 

4 TEM3EA, Tf MSB, TEM2Bd, TEM2AE ,F ,TcK2liA, TEM4BA, uN jM» 

R FCRMATd • ,G15.6 , l£X,G15.6 ) 

10 FCRMATd ‘ ,4G15.6 ) 

CC 960 J=1,NSTGT 

3RE = DLUG(CSQPT(CUR( J)*-*2»CR I( J »>K*2J )/T:.AMP 
SIM = CATAMCWI ( J)/CW?.{ Jl )/TSAMP 
IF <CM(J) .EQ. 0.00) GC TC 4974 
WRITE(6,1C) CWR{ J) ,CWKJ) ,SRE,SIM 
GC TO 960 

4974 WRITEIc^E) CWR(J),SRE 
960 CCNTINUE 
GC TG 1 
A6 RETURN 
END 
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c 

c 

, SLBBOUTINt CNTRL (KPLAM,NniG,NC,Ni:T,.^S,NS7QT,NSrjT2» 

1 ICWT , INT ,G.J,U9, C^lR^C^^ I , T , AC L ,KVEC i VcCRNt ViCIr4,VECR» 

2 vrc I ,VRV,VlV,V^.RV,VRIV,AFPI'•1E,bPKiM!£,P,0,Cl,C/l,C3,MltM2fM3,^»4, 

3 M5fNic,A,e,C,L),TE«5HA,TEMA;,TEM2AA, T=MdA,IcMAtj,Vi,\/2t 

4 TEM3EA,TEMBB,TEf^2BB, TtR2AH ,Ff T£t-l<.iiA, riiM4BA,Ui-40M) 

IMPLICIT REAL (K,Ltl^) 

REAL«E F, ICNT, INT ,GN,CS ,ACL 
REALMS SR:,SIM,CUI ,CWR 
R£AL*8 OLCG,DSOR7,DATAN 
CCMMGN/T/TSAMP,CFCLES,CZEPCS 

01 MEN SIGN ICNT(N5T0T2 ) , INT ( NS TuT2 > , GN (nSTO f, rtS TjT J , 09 (N ST0T2 ) f 

1 CWR (NSTOT) ,CWI (NSTQT ) , T I N V (NDIG ,.Jii I G » t ACL IN j T uf , F JT » , 

2 RVcC (NSTCTflFSTOT » , Ve CRN ( N STOT ,NS F JH ,VtCiN(KiF JF, .<STUT ) , 

3 VtCR (NSTCT.NSTGT ) , VECI (NSTGF.NSTjF ) , V-.V ( NSTOF ) , VI V ( NiTOT ) , 

4 VRKV(NSTCT) , VR I V ( MSTC T ) , A F R I MciNS , NS ) , QPkI Nt ( hS tNC ) . 

5 P( NS, NS) ,G«NS ,NC) ,C1 (NC.NS F.CZlNL.uDIG) ,C3(NC ,NS) 

OIMENSICN MKNOIG.I'ICIG) , M 2 ( NDI GiNOi G ) ,M3 (GDI G, i«i;> ) , M4 ( NOI G,MS) , 

1 M 5 ( rj C I G ,N S ) , M6 ( N C I G , N S ) , A ( NPL AnF , NPlAi;T ) , ii< NP LAUF , NC ) , 

2 C( NACT.NACT) ,0(NACT ,NC ) , Tr^/5RA(^^0I G,;4S » , FEHaa 1.'IS,NS) , 

3 TEM2AA(NS,NS) ,FEMBA(NCIG,NS),FcMAi(NS,NDIi,) ,VA(USFJF2) , 

4 V2(NST0 72),TEM3eA(NDIG,NS) ,f E MB8 ( .,0IG, tjul G) , F cl-i2BtJ ( iv 0 I G, MO IG ) , 

5 TEH2A6(NS,N0IG) , F ( NSTOl , N570T ) ,TtM2GA( NOlG.nS ) , TtMAb A ( MDIC ,NS ) 

6 .DNaR^F ACT.NC) 

C 

CALL CIFFEC{NPLANT,NOIG,NC,NACr,N5,NSFGT,N^FuT2, 

1 ICNT , INF ,GN,09, C^R,CKI,TINV,ACL,RVLC,VcCRi>l, VbCIi<, VECR, 

2 VeCI ,VRV,vrV,VR!7V,VRlV,AFKIME,BP\lME,P,G,Cl,k,L,C3,M.l ,M2,M3,M4, 

3 M5,M6,A,e,C,D,TEM56A,TEMAA,TEM^;AA, FcM3A,FtMAtj,Vl,V2, 

4 TFM3eA,7EMB8,TeM2e9,TeN2Ae,F,TtM2aA,FcM4SA,0NuM) 

C 

CO 91 1= l.NPLANT 
DC 90 J= I.NPLANT 

90 APR I,ME( I ,J)= A(I,J) 

CC 91 J= 1,NC 

91 APR IMEn,j4-NPLAN7)= E(I,J) 

CO 94 1= 1,NACT 

CC 94 J=1,NACT 

94 APPIM£( I+NPLAMT, J*NPLAN7)= C(I,J) 

C 

DC 96 I=1,NACT 
00 96 J=1,NC 

96 BPK IME ( I iNPLAfn, J)=0( I , J) 

CALL FGPC (APRIM6 ,BPBIFE,P ,C,NS ,NC,7SAMP) 

CALI. MMUL7(C3,AP«IME,7FR5EA,MC,NS,iJ;>) 

CALL NACC (7lM5BA,C1,7E^'5BA,KC,NS,uI 
CALL NMULT(C,7EM5eA,7ef'AA,NS,NC,iFS) 

C.^U. VACT (TF/-1A A ,P ,7F«AA ,NS ;NS ,0! 

CALL f'MLL7(0,C2,7EMAB,N5,NC,NDIG) 

CALL ^A7I^( 7I,4V,NCIG, 1 ) 

. CALL MACCI7 Ii'IV,M2,TI')V,N0IG,NDIG. 1) 

CALL MIMV (MDIG*MO!G,T INV, m: IG ,01, VI , V4l) 

CALL MMUL7I M6, A^R !ME ,T£,*'2P.A ,MD IG,;io ,fj>) 

CALL NACC(7EiV2BA,F4,TfcP2BA,M0IG,Mi,0» 

CALL MVUl,7( ^EV2PA,TeMAA,rE^3BA,l<Di^F,:-^S,l'^S) 

CALL MMOL r (M5, APR !M;5 , T EMAE A ,iiD IG,i;s ,.'IS) 

CALL V A D n ( T E M 4 P A , T “ 1^3 e A , T E (« 3 B A , N I » i G , M S , 0 ) 

CALL MArD<7'iM33A,M3,T£>'’';A,i\nlG,N.Mv>» 

CALL M.-’ULliriHVfl'.Vlit:/ , ; ,f 

CALL MMUL7 (Q,C2, TE'^ZAE , NS,, VCrrulo) 
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CALL MMUL7(TEM2BA,TEM2AB,TLf'28ri,.NJiG,NS,Ni)iG) 

CALL yACC{TcM2P3,Ml,rEf'29e,f<0IG,.NJIG,U) 

CALL NPULT(TINViTEM2eefTEPEB,NUlG,(MUlGfNGlG) 

DC 125 1=1, NS 
CO 125 J=1,NS 
i25 F( I , J)=TcMAA( I , J» 

CG 76 I=l,NOIG 
CO 76 J=1,NS 
F( J, I +NS ) = TfPA8 (J , I ) 

76 F( UNS, J)=TEMBA{ I , J» 

DO 77 I=1,NDIG 

Ce 77 J=1,NDIG 

77 F( H-NS, J + 6S)=TEMRR(I ,J) 

C** 

CALL EALAAC ( N STC T , NS TCT , F , ILOW , I N I GH , D9J 
CALL ELPFES ( NSTCT , NSl CT , I LCW , I HI GH , F , I NT ) 

CALL hCR2 ( NSTQT , N STU T , I LCU , I HI GH, F , CWa , C W l , AC L, iCNT , C46) 

A6 rtFTURN 
END 


SUB BOLT I NE LEAD! XB , Y I ,NPCL£S, POLtjNZcKJS , ZcBU, GAIN.CJMPP.RHT , 
1 CCMPZ.SF) 

CIPENSICN POLE (2»6) ,ZFRC(5,8I 
ANGLE=0. 

GAIN=1. 

IF( NPCLES .EU.O) GO TO 15 
CC 10 I=1,NP0LES 

GAIN=GAIiN>»SUKT( ( XF-PULE ( 1 , I ) ) *=»2+ ( Y I -POLt U,iJ >**Z) 

TcMP = ATAFil(YI-P0LE(2, I) )/{XP-PULtU,I») ) 

IFIPOLS ( Ir D.GT. XR) T EPP= 3 . 14 1592<-T UMP 
10 ANGlE=ANGLL*-TEMP 
16 F:rPAT(7F10.2J 
15 IF(NZEBCS.EO.O) GC to 25 
CO 20 Irl.NZfcROS 

GAIN=GAIF/«S0RT( (XR-Z£PC( 1, I) YI-ZcR0(2,I n»*2)J 

T£MP=ATA6( (Yl-ZERr(2 , I) I / ( XR-Z=RjI i , I ) ) » 

IF(ZEFQ(1,I).GT.XH) T EF F = 3 . 141 592^-T lMP 
20 A.NGLt = A.;GLF-T£MP 

25 PHI = ANGLr«-3. 141592 

FHiZ = ATANlYl/(XR-CCMPZn 
IFICCMPZ.GT .XR) Ph! Z= 3 . 14 1 592 ^►PHl Z 
CCMFP = XP. - YI/TAN( FF I Z-FF I ) 

GAIN=GAlF^ SOP.TI (XR-CCFPP)«'*24-(YI)*'«‘2t 
GA IN = SF*GAI,N/ ( SCRT( (XP-CCMPZJ-^OZMYD^+'Z) ) 
GAIN^-GAIN'^d.-CCMPZ)/ (i.-CCPPP) 

RETURN 

END 
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SLBPGITINE DMAT( C,D^C^ t ^ACT f.NC,OAIN I 
DI^cMSIa^ DNGM(NACT ,NC) ,0(NACT,NC) 

CC 10 I=ltNACT 

DO 10 J=1,NC 

C( I , J )=DNCM( I ,J )«GAIN 

anitRN 

END 


SOBROt.tTINE SIZ5(KPLAN'r,Nr'IG,!:C,Wc'!’) 

CCCCC 

CCCLX SIZE IS A !JSFR-VR rT?r?: PPOnmt> DE‘"'EP^*T*rFS THE SIZE OP 

CCCCC THF SYSTF*'. I? SHOlIin PE PROGRAEKEn AE IN '"HTS EXAMPLE. 

CCCCC MPLANT = T’!'- NHMPtE OF DIFFERENTIAL EQUATIONS OP TH^' PLAHT 
CCCCC ?.'DTG = THE N'H.'inFR r>^ DTFFEPENCE EOUATIOH TH THS CONTROL 
CCCCC MC = THE NHMBF” OF CONTROLS 

CCCCC NACT = THE NHMDEP OF ACTHATOR PI FEE RFN'^T A L EgUATTOES. 

CCCCC NOTE: THIS PEAKS THAT NPLANT DOES NOT INCLUDE ANY 

CCCCC ACTUATOR DYNAMICS 

CCCCC THE PROGRAM SHOWN HERE IS FO? THE LOMGITnniTUL ATTITUDE 
CCCCC CONTROL-KHELL STEERING PROBLFP DISCUSSED IN THE BODY OF THE 
CCCCC REPORT 
CCCCC 

NPLANT=3 

NCIG=3 

HACT=1 

NC=1 

RETURN 

END 

C 

C 

c 

c 

c 

c 

c 

c 

c 

c 

S OP.F.OUTINE PLANT (JJPLn NT, NDIG,NC, NACT, US , ESTOT ,N STO T2 , 

1 TCNT,1N'’, GN, DP,CWR,CWT,TINV , ACL , R VEC , YEC P»» , VECTN , V^CP , 

2 VFCI,VPV, VTV,V^'V,yriV, APRIME,PPRTMF,P,0,C1,C2,C3,Ml,M2,H3,nU, 

3 «5,n6,A,B,r,D,TEMSBA,TEMAA,"’EM2AA,TEM'3A,TFMAB,V1,V2, 

« T EM 3 E A , T Ef* BB , T EE 2BB , T FH 2A n , F , TE M2U A, '^EM« P A, DNO M) 

IMPLICIT REAL (K,L,M) 

PEAL*F. ?, ICNT,INT,GN,D9, ACL 
R’=’AL*H EF F,SrM,CWr ,CWR 
REALf^a DLOG,DS0PT,DATAN 
COM«ON/T/TS AMP, CP'^LFE , C7.EPOS 
OIMEN3TON Z’^F 0 {2 , R ) , POLE (2 , F) 

D.Tt'ENEIOH ICK'’’ (NSTO'^2) ,INT (NSTOT2) ,GK (NST0'^,MST0’') ,nU fN3TOT2) , 

1 CW'> (KSTOT) ,C WI (NSTOT) T N V r ND IG , MD TG) , ACL f NSTOT, K S'T’OT) , 

2 R VEC (NSTOT,KST'^T) , VECRN (NSTOT, KSTOT) , VECIN {NSTOT, KSTOT) , 

3 VECP (KS”'0'",NSTO™) , VEC I (N S"'" KSTO'") , VRV (HS "OT) ,VIVfNSTOT), 

U VP? V (KSTO’') , V^I V (KSTOT) , A PPI MF ( KS , NS) , BP F I M E ( N S , N C) , 

5 P(?;S,NS) ,Q(fr3,NC) ,C1 (VC, ’IS) ,C2(KC, "DIG) ,C3 C’C,NS) 

DTM’^NSTON 1 (KD'G, NDTr.) , M 2 ( VDTG , KD’G) , F 3 (’JDI G , N S) , MU <’NDTG,NS), 

1 M5(?:DTG,M?) (MDTG,VS) ,A (FPIAHI, KFLAVT) , BfNPL AK"' ,»’C)., 

2 C(HACT,NAvC'") , n (VAC^, VC) , " EB A ( VD ’'G , NS) , v A A (IIS , VS ) , 

3 TFr2 AA (NS,HS) ,'PFMRA (MD'G, ?!F) ,TEMAB (NS,"DIG) , VI (M3TOI2) , 

4 V2 (HS"T)T?) ,TE:MPA (I;DTG,N3) ,TFMPP(NDIG,?’DTG) (VDIGjMBIG) , 

' TFfi2AP(VS, VDIG) , F(VS"^0'", rxTC'") , "F M 2R A (NPIG ,KS) ,"EFU BA f NDIG , NS) 

6 ,DKOM (Nr,c?,wr) 

CCCCC 

CCCCC i' T A MT IS A HS'^F -- W~ CTT'^’- PROGFAM THAT DETERMIKEG THF 
CCCCC DIFFEREVTTM. FvT'ATI"'! '’.'■F IRH PL7,NT, MOD" P *^01 F IC A T f Y IT 
CCCCC DEFl'-FS THF. MA'rICF:' .ALP DNO«. TH F ^ IH 1 1 V F C.ALI.ED, . 

CCCCC Oi;CF -PFF .Hi.i'I, AND T'!*' r F SO .'•'AY CC-VTAT'? READ ST ) ~E'"r'F'’'S .AND 
CCCCC THE LIKE. 
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ccrcc THF P^onPUM r,F!OW” TS EOF TEE LONGT T!1PIN AL ALTITUDE 

CCCCC CCNTROI-WIIELL STEFPTNG PHOBIFF' n.TSCUSSfD IN THE non? OF THE 

ccccc p.EPOirr 

CCCCC 

A (1, 1) = -6. 2a2a7.E-1 
A (2,1) =-?.27919E-3 
A (1,2) =207.064 
A (2,2) =-4 . uqPSUE-l 
A(A,2)=1. 

B (1 , 1) =-1 . 408fl0F-7 
B (2, 1) =-1 . 7392aE-2 
C(1,1)=-10, 

DNO"! (1,1) =10. 

RETUPH ■ 

EN P 
C 
C 
C 
C 

C - 

c 

c 

c 

c 

c 

S UBEOtlTIN^ DIFF^O ( NPLA NT , I'DTG , NC , K ACT , M5 , NSTOT , NSTOT2 , 

1 ICMT, TNT, G^’, P«,CVR ,CWI, TTNV , AC L ,Evrc , VECPN , V FC IN , V’^CP , 

2 VECT, V?V, VI7,VF-’7 , VPIV, A PPI ME, B PF^ !1F ,F ,0 ,C 1 , C2 , C3 , M 1 , M2 , M3 , N4 , 

3 M5,H6,S ,B,C,D,'^rM53A,T^MAA,TE«2AA,TFrBA,TSMAB, V1, V2, 

4 TE"3 BA,"- E''PP,':FM 2BE ,TE.M 2AE , ? ,TFM2BA, •TFM4BA, DNOM) 

IMPI.ICIT FEAL (F,L,K) 

PEAI.^9 F, ICNT ,INT,GN,D9, ACL 
??AI«8 SPE,FlM,CyT,CHR 
REAI.’*Q DLOG, PSOPT, DATA?! 

CO F«ON/T/"'S A M P, CPO L ES ,C7ER0? 

DTK^NSrOH TCNT(NSTOT2) ,rN?(N?TOT2) , GK NSTOT) , P9 (K STOT2 ) , 

1 CWP (FFTOT) ,C'WIC* SrOT) , TIE V ( vpxG , NDTG) , ACl, (TETOT ,NSTOT) , 

2 RVEC (N??O'^,NST0?) ,VFCRN (ESTOT, KS’^OT) ,vrciN (K SCOT , NSTOT) , 

3 VECF (N3TOT ,NSTOT) ,V^CT (J:FT0T,’:STDT) , V- V (HSTO'^) ,VTV(KSTOT), 

4 VPPV (KSTO-") , VRI7 (N’STOT) ,APFIFE (I:E,KF) ,PPRIMF (HE ,KC) , 

E r {”3, NS) ,0 (EE,FC) ,C1 (!IC, KE) ,C2 (NC, HPIG) ,C3 (NC, NE) 

DIF f; ME ION F1 (MD"G, NPIG) , .M2 (KDIG, NDIG) ,E? (KDIG,ME) , ,M4 (NPIG,N.S) , 

1 M.E (!J DIG, ?:.3) , F6 (’;r iG,MP) ,A (FPIA NT, FFI A*’T) ,B (FPL AN*" ,r;C) , 

2 C (FACT, MAC") ,D (N .ACT, NC) ," F M 5B A ( NPTG , NS) , TEE A A (ME , N.E ) , 

3 TFM2 A.A(VS,NE) , TE MBA ( N P"G , fIS ) ,TFMAB (?.’S,VDIG) ,V1 (MSTOT2) , 

4 V2(v?'TOT 2) ,'T'EM3PA (M PIG , M S) , TE MB B (MOTG, NDTG) , TF.^ 7 pp ( Mpir, , noIG) , 

5 TFM2.A" (KS,NDrG) , (H STO T , NSTOT) , TF K 2B A ( MDIG , HE) ,T EM4BA (MDIG, NE) 

6 , D.MOM (NACT, 'IC) 

CCCCC 

CCCCC PIFFEO IS A HE EB - WP ITTEK PRCGFAM TBA” DFTFRMIFFE T^F DTFFEPFNCF 
CCCCC ='QL'ATtrNS USED IN THE DIGITAL CONTRCL. «OPF EP EC I F ICTAL LY 
CCCCC I F DEFINES THE MATRICES: M 1 , M2, M .3, M 4 , ME , M6, AMD C1,C2,C3 

CCCCC IN "EPMS OF "SAMP AND ANY CONTPOT C0HN"ANTS, TT IS CALLED 
CCCCC SEVERAL TIMFS IN GETTING TH'^ COKrENEA"'OM (rVS = 2) *^0 CONTFOL 

CCCCC P AFA.M FTEF.E SM07LD NOT BE READ IN Tfl.E .ENP^ONIJNE , B'lT r'fnVLD 
CCCCC INST'^AP BE KEAD AT ""E DEGI’’rT’TG OF Til" PPOGPAM .AMP PASSED 
CCCCC TO FIFF'^O TH .3 I.f- 1”^ i,:!: D COM.'iOM. T!'^ ECNAT"OK OF "BE 
CCCCC CO''D"':?ATT0M IC'ST priGnPA^HED HESF, AMD "HNS, NDTG, "II?. 

CCCCC M’lNpr.B OF Di ? p CFC E FOIIA"IO!IS M US " T’lCLUDF IRb CO M? FNS A TTON 
CCCCC coin "lev, 

CCCCC "HE hcogbaM SHOVN '/.ere is FOR THE I ON GIT'/ PIN AI. ATTITflDE 
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CCCCC CO«TPOL-WHi;i,I, STKFRI'Jr, PF03LFM DISCUSSED IN THF BODY OF THE 
CCCCC PFPORT 

ccccc 

F01F = PXP (-TSS iMP/16. ) 

E016 F3=S7. 2'-''5 77«S1'«' (1 . f F016)/2, 

C2(1,3)=1. 

C Wf, ‘5110IIT 

« 1 (1 , 1) = E016 
?*3(1,2) =-FO-!6F3*2. 16 
NU (1 , 2 ) =E016F3*2. 16 
C COMPEMSMION 

TT2=k:XP (-?EANP*CPOLSS) 

=''XP (-TS7»riP*CZrROS) 

M1(2,2)=TT2 

r.LL= M.-'^T2) / n.-'^'f’l) 

N 1 (2, 1) =-rl!,*TT1 
N2(2,1)=GLL 

M 3 (2, 3) --r,LL*TT1*247. 5 3 
N4 (2, 3) =GLL*247. 53 
C DELAY 

wl(3,2)=1. 

FETOP.H 

END 
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The utility routines are: 


MMULT 

MADD 

MATIN 

MATEQ 

SCAMAT 

FGPQ 

MINV 

CDIV 

HQR2 

BALANC 

ELMHES 


matrix multiplication 

matrix addition 

matrix initialization (0 or I) 

equate two matrices 

matrix-scaler multiplication 

calculates <p given F 

matrix inverse 

complex nvimber division 

eigenvalue solver (QR algorithm) 

preprocesses matrix for HQR2 

postprocesses matrix from HQR2 
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